




*********************************Bootstrap Table 1********************************************************************
*********************************Bootstrap Table 1********************************************************************
*********************************Bootstrap Table 1********************************************************************

*-----------------------------------------------------------------------------------------------------------------------
* Part I: The following code is to replicate the "bootstrap confidence interval" in Table1 (columns 1, 2 , 5, 6 )
*-----------------------------------------------------------------------------------------------------------------------




 clear mata
   mata:

void GMM_DLW_CD(todo,betas,crit,g,H)
  {
	PHI=st_data(.,("phi"))
	PHI_LAG=st_data(.,("phi_lag"))
	Z=st_data(.,("const","m_lag","l_lag","y1","y2","y3","k","k_lag","exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))
	X=st_data(.,("const","m","l","x1","x2","x3","k","exp_dummy99","exp_dummy00","exp_dummy01","exp_dummy02","exp_dummy03","exp_dummy04","exp_dummy05"))
	X_lag=st_data(.,("const","m_lag","l_lag","x1_lag","x2_lag","x3_lag","k_lag","exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))
	Y=st_data(.,("q"))
	C=st_data(.,("const"))
	

	EXPORT=st_data(.,("exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))

	OMEGA=PHI-X*betas'
	OMEGA_lag=PHI_LAG-X_lag*betas'
	OMEGA_lag_pol=(C,EXPORT,OMEGA_lag)
	g_b = invsym(OMEGA_lag_pol'OMEGA_lag_pol)*OMEGA_lag_pol'OMEGA 
	XI=OMEGA-OMEGA_lag_pol*g_b /*XI=epsilon: recover epsilon_it in equation (23): shock in productivity;*/
	crit=(Z'XI)'(Z'XI)
  }
///
///
void DLW_TRANSLOG()
	{
	initialvalue=st_data(1,("initialConst","initialm","initiall","initialx1","initialx2","initialx3","initialk","initialexp_dummy99","initialexp_dummy00","initialexp_dummy01","initialexp_dummy02","initialexp_dummy03","initialexp_dummy04","initialexp_dummy05"))
	S=optimize_init()
	optimize_init_evaluator(S, &GMM_DLW_CD())
	optimize_init_evaluatortype(S,"d0")
	optimize_init_technique(S, "nm")
	optimize_init_nmsimplexdeltas(S, 0.1)
	optimize_init_which(S,"min")
	optimize_init_params(S,initialvalue)
	p=optimize(S)
	p
	st_matrix("beta_dlwtranslog",p)
	}
end

capture program drop dlw_translog
program dlw_translog, rclass
preserve
**********************************************
sort firm year
mata DLW_TRANSLOG()
end 
use "D:\Nanjing\2019\pollution2\submission\data\boot1",clear
rename frdm firm
gen hangye=int(cic_adj/100)
drop if hangye<13
drop if hangye>42

*keep if hangye==27|hangye==14
capture drop cic2
capture drop cic4 
gen cic2=floor(cic_adj/100)
egen gcic4=group(cic_adj)
egen SDcic=sd(cic2), by(firm)
egen cic_y=group(cic year)
gen Productivity=tfp_acf_new
*----------------------------------------------------------------------
* bootstrap starts from here
*----------------------------------------------------------------------

capture program drop prodboot
program define prodboot, rclass
preserve






* Re-estimate firm-level production function for each sample from bootstrap
do "D:\Nanjing\2019\pollution2\do file\for accre\do file\tfp_combine.do"

xi:heckman ratio1 dyear*, noconstant twostep sel(exp_dummy=  tfp_acf_new clean_dummy_so2_lag clean_dummy_dust_lag clean_dummy_cod_lag dcic2* lnv_lag lnkl_lag dprov* owner* dyear*)  

gen coef1_bs=_b[dyear1]
gen coef2_bs=_b[dyear2]
gen coef3_bs=_b[dyear3]
gen coef4_bs=_b[dyear4]
gen coef5_bs=_b[dyear5]
gen coef6_bs=_b[dyear6]
gen coef7_bs=_b[dyear7]
capture drop exp_dummy2
gen exp_dummy2=.
replace exp_dummy2=exp_dummy*coef1_bs*dyear1 if year==1999
replace exp_dummy2=exp_dummy*coef2_bs*dyear2 if year==2000
replace exp_dummy2=exp_dummy*coef3_bs*dyear3 if year==2001
replace exp_dummy2=exp_dummy*coef4_bs*dyear4 if year==2002
replace exp_dummy2=exp_dummy*coef5_bs*dyear5 if year==2003
replace exp_dummy2=exp_dummy*coef6_bs*dyear6 if year==2004
replace exp_dummy2=exp_dummy*coef7_bs*dyear7 if year==2005
drop  coef1_bs coef2_bs coef3_bs coef4_bs coef5_bs coef6_bs coef7_bs






* For so2 + OLS: column 1 of Table 1
areg lnso2  exp_dummy2 so2_new1 Productivity lntariff ///
 i.year i.prov i.ownerships, a(cic) cluster(cic_y)

return scalar exp_t2c1=_b[exp_dummy2]
return scalar abt_t2c1=_b[so2_new1]
return scalar tfp_t2c1=_b[Productivity]
return scalar tar_t2c1=_b[lntariff]

*************************************************************************************************
* For so2+ IV: column 2 of Table 1
areg lnso2  exp_dummy2 so2_new1  ///
Productivity lntariff lnkl_lag lnso2_ini  i.prov i.year,a(cic) cluster (prov_c)

return scalar exp_t2c2=_b[exp_dummy2]
return scalar abt_t2c2=_b[so2_new1]
return scalar tfp_t2c2=_b[Productivity]
return scalar tar_t2c2=_b[lntariff]
return scalar tar_t2c2=_b[lntariff]
return scalar klr_t2c2=_b[lnkl_lag]
return scalar iei_t2c2=_b[lnso2_ini]

********************************************************************************************************

* For dust+OLS: column 3 of Table 1
areg lndust  exp_dummy2 dust_new1 Productivity lntariff /// 
  i.year i.prov i.ownerships, a(cic) cluster(cic_y)

return scalar exp_t2c3=_b[exp_dummy2]
return scalar abt_t2c3=_b[dust_new1]
return scalar tfp_t2c3=_b[Productivity]
return scalar tar_t2c3=_b[lntariff]

************************************************************************************************************************************
* For dust +IV: column 4 of Table 1
areg lndust  exp_dummy2 dust_new1 ///
  Productivity lntariff lnkl_lag lndust_ini  i.prov i.year, a(cic) cluster (prov_c)

return scalar exp_t2c4=_b[exp_dummy2]
return scalar abt_t2c4=_b[dust_new1]
return scalar tfp_t2c4=_b[Productivity]
return scalar tar_t2c4=_b[lntariff]
return scalar klr_t2c4=_b[lnkl_lag]
return scalar iei_t2c4=_b[lndust_ini]

***************************************************************************************************************************************


end

capture drop indc
egen indc=group(firm)
drop if indc==.

**************************************************************************************************
* The following code generate columns 1, 2 , 5, 6 in Table 1.
*************************************************************************************************

bootstrap  exp_t2c1=r(exp_t2c1) abt_t2c1=r(abt_t2c1) tfp_t2c1=r(tfp_t2c1) tar_t2c1=r(tar_t2c1)  ///
 exp_t2c2=r(exp_t2c2) abt_t2c2=r(abt_t2c2) tfp_t2c2=r(tfp_t2c2) tar_t2c2=r(tar_t2c2) klr_t2c2=r(klr_t2c2) iei_t2c2=r(iei_t2c2)   ///
 exp_t2c3=r(exp_t2c3) abt_t2c3=r(abt_t2c3) tfp_t2c3=r(tfp_t2c3) tar_t2c3=r(tar_t2c3)   ///
 exp_t2c4=r(exp_t2c4) abt_t2c4=r(abt_t2c4) tfp_t2c4=r(tfp_t2c4) tar_t2c4=r(tar_t2c4) klr_t2c4=r(klr_t2c4) iei_t2c4=r(iei_t2c4)   ,  ///
 reps(50)  seed(1) nodrop cluster(indc) idcluster(newid): prodboot
 


 
 
 
 
 
 
 
 

*-----------------------------------------------------------------------------------------------------------------------
* Part II: The following code is to replicate the "bootstrap confidence interval" in Table1  Column 3 - 4 and 7-8
*-----------------------------------------------------------------------------------------------------------------------
   clear mata
   mata:

void GMM_DLW_CD(todo,betas,crit,g,H)
  {
	PHI=st_data(.,("phi"))
	PHI_LAG=st_data(.,("phi_lag"))
	Z=st_data(.,("const","m_lag","l_lag","y1","y2","y3","k","k_lag","exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))
	X=st_data(.,("const","m","l","x1","x2","x3","k","exp_dummy99","exp_dummy00","exp_dummy01","exp_dummy02","exp_dummy03","exp_dummy04","exp_dummy05"))
	X_lag=st_data(.,("const","m_lag","l_lag","x1_lag","x2_lag","x3_lag","k_lag","exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))
	Y=st_data(.,("q"))
	C=st_data(.,("const"))
	

	EXPORT=st_data(.,("exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))

	OMEGA=PHI-X*betas'
	OMEGA_lag=PHI_LAG-X_lag*betas'
	OMEGA_lag_pol=(C,EXPORT,OMEGA_lag)
	g_b = invsym(OMEGA_lag_pol'OMEGA_lag_pol)*OMEGA_lag_pol'OMEGA 
	XI=OMEGA-OMEGA_lag_pol*g_b /*XI=epsilon: recover epsilon_it in equation (23): shock in productivity;*/
	crit=(Z'XI)'(Z'XI)
  }
///
///
void DLW_TRANSLOG()
	{
	initialvalue=st_data(1,("initialConst","initialm","initiall","initialx1","initialx2","initialx3","initialk","initialexp_dummy99","initialexp_dummy00","initialexp_dummy01","initialexp_dummy02","initialexp_dummy03","initialexp_dummy04","initialexp_dummy05"))
	S=optimize_init()
	optimize_init_evaluator(S, &GMM_DLW_CD())
	optimize_init_evaluatortype(S,"d0")
	optimize_init_technique(S, "nm")
	optimize_init_nmsimplexdeltas(S, 0.1)
	optimize_init_which(S,"min")
	optimize_init_params(S,initialvalue)
	p=optimize(S)
	p
	st_matrix("beta_dlwtranslog",p)
	}
end

capture program drop dlw_translog
program dlw_translog, rclass
preserve
**********************************************
sort firm year
mata DLW_TRANSLOG()
end 
use "D:\Nanjing\2019\pollution2\submission\data\boot_app2",clear
*
*----------------------------------------------------------------------
* bootstrap starts from here
*----------------------------------------------------------------------
*keep if _merge==2|_merge==3
capture program drop prodboot
program define prodboot, rclass
preserve






* Re-estimate firm-level production function for each sample from bootstrap



xi:heckman ratio1 dyear*, noconstant twostep sel(exp_dummy=  tfp_acf_new clean_dummy_so2_lag clean_dummy_dust_lag clean_dummy_cod_lag dcic2* lnv_lag lnkl_lag dprov* owner* dyear*)  

gen coef1_bs=_b[dyear1]
gen coef2_bs=_b[dyear2]
gen coef3_bs=_b[dyear3]
gen coef4_bs=_b[dyear4]
gen coef5_bs=_b[dyear5]
gen coef6_bs=_b[dyear6]
gen coef7_bs=_b[dyear7]
capture drop exp_dummy2
gen exp_dummy2=.
replace exp_dummy2=exp_dummy*coef1_bs*dyear1 if year==1999
replace exp_dummy2=exp_dummy*coef2_bs*dyear2 if year==2000
replace exp_dummy2=exp_dummy*coef3_bs*dyear3 if year==2001
replace exp_dummy2=exp_dummy*coef4_bs*dyear4 if year==2002
replace exp_dummy2=exp_dummy*coef5_bs*dyear5 if year==2003
replace exp_dummy2=exp_dummy*coef6_bs*dyear6 if year==2004
replace exp_dummy2=exp_dummy*coef7_bs*dyear7 if year==2005
drop  coef1_bs coef2_bs coef3_bs coef4_bs coef5_bs coef6_bs coef7_bs







*************************************************************************************************
* For so2+ OL: column 3 of Table 1

 ivreghdfe lnso2  Export Abatement_so2 Productivity Tariff lnKL , a(indc_firm year cic)

return scalar exp_t1c2=_b[Export]
return scalar abt_t1c2=_b[Abatement_so2]
return scalar tfp_t1c2=_b[Productivity]
return scalar tar_t1c2=_b[Tariff]
return scalar klr_t1c2=_b[lnKL]

* Column  7 of Table 1
capture drop ub
egen ub=pctile(lndust), p(97)
 ivreghdfe lndust  Export Abatement_dust Productivity Tariff lnKL if lndust<ub, a(indc_firm year cic)

return scalar exp_t2c2=_b[Export]
return scalar abt_t2c2=_b[Abatement_dust]
return scalar tfp_t2c2=_b[Productivity]
return scalar tar_t2c2=_b[Tariff]
return scalar klr_t2c2=_b[lnKL]
********************************************************************************************************


*----------------------------------------------------------------------------------------------------
* Column 4 of Table1
ivreghdfe lnso2_d  exp_dummy2 so2_new1_d  ///
 tfp_acf_new_d lntariff_d lnkl_d, a(cic prov year ownership)

return scalar exp_t3c4=_b[exp_dummy2]
return scalar abt_t3c4=_b[so2_new1_d]
return scalar tfp_t3c4=_b[tfp_acf_new_d]
return scalar tar_t3c4=_b[lntariff_d]
return scalar klr_t3c4=_b[lnkl_d]

*----------------------------------------------------------------------------------------------------
* Column 8 of Column 1
ivreghdfe lndust_d  exp_dummy2 dust_new1_d  ///
 tfp_acf_new_d lntariff_d lnkl_d , a(cic prov year ownership) 
return scalar exp_t4c4=_b[exp_dummy2]
return scalar abt_t4c4=_b[dust_new1_d]
return scalar tfp_t4c4=_b[tfp_acf_new_d]
return scalar tar_t4c4=_b[lntariff_d]
return scalar klr_t4c4=_b[lnkl_d]


end

capture drop indc2
egen indc2=group(firm year)
drop if indc2==.


**************************************************************************************************
* The following code generate columns 3, 4 ,7, 8 in Table 1.
*************************************************************************************************
bootstrap   exp_t1c2=r(exp_t1c2) abt_t1c2=r(abt_t1c2) tfp_t1c2=r(tfp_t1c2) tar_t1c2=r(tar_t1c2)  klr_t1c2=r(klr_t1c2) ///
exp_t2c2=r(exp_t2c2) abt_t2c2=r(abt_t2c2) tfp_t2c2=r(tfp_t2c2) tar_t2c2=r(tar_t2c2)  klr_t2c2=r(klr_t2c2) ///
exp_t3c4=r(exp_t3c4) abt_t3c2=r(abt_t3c4) tfp_t3c2=r(tfp_t3c4) tar_t3c2=r(tar_t3c4)  klr_t3c2=r(klr_t3c4) ///
exp_t4c4=r(exp_t4c4) abt_t4c4=r(abt_t4c4) tfp_t4c4=r(tfp_t4c4) tar_t4c4=r(tar_t4c4)  klr_t4c4=r(klr_t4c4), ///
 rep(50)  seed(5) nodrop : prodboot





************************************






 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
*********************************Bootstrap Table 3********************************************************************
*********************************Bootstrap Table 3********************************************************************
*********************************Bootstrap Table 3********************************************************************

 
 



*-------------------------------------------------------------------------------------------------
*Part I: The following code generate "bootstrap confidence interval" in Table 3 column 1, 2, 5, 6,
*------------------------------------------------------------------------------------------------
 

 clear mata
   mata:

void GMM_DLW_CD(todo,betas,crit,g,H)
  {
	PHI=st_data(.,("phi"))
	PHI_LAG=st_data(.,("phi_lag"))
	Z=st_data(.,("const","m_lag","l_lag","y1","y2","y3","k","k_lag","exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))
	X=st_data(.,("const","m","l","x1","x2","x3","k","exp_dummy99","exp_dummy00","exp_dummy01","exp_dummy02","exp_dummy03","exp_dummy04","exp_dummy05"))
	X_lag=st_data(.,("const","m_lag","l_lag","x1_lag","x2_lag","x3_lag","k_lag","exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))
	Y=st_data(.,("q"))
	C=st_data(.,("const"))
	

	EXPORT=st_data(.,("exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))

	OMEGA=PHI-X*betas'
	OMEGA_lag=PHI_LAG-X_lag*betas'
	OMEGA_lag_pol=(C,EXPORT,OMEGA_lag)
	g_b = invsym(OMEGA_lag_pol'OMEGA_lag_pol)*OMEGA_lag_pol'OMEGA 
	XI=OMEGA-OMEGA_lag_pol*g_b /*XI=epsilon: recover epsilon_it in equation (23): shock in productivity;*/
	crit=(Z'XI)'(Z'XI)
  }
///
///
void DLW_TRANSLOG()
	{
	initialvalue=st_data(1,("initialConst","initialm","initiall","initialx1","initialx2","initialx3","initialk","initialexp_dummy99","initialexp_dummy00","initialexp_dummy01","initialexp_dummy02","initialexp_dummy03","initialexp_dummy04","initialexp_dummy05"))
	S=optimize_init()
	optimize_init_evaluator(S, &GMM_DLW_CD())
	optimize_init_evaluatortype(S,"d0")
	optimize_init_technique(S, "nm")
	optimize_init_nmsimplexdeltas(S, 0.1)
	optimize_init_which(S,"min")
	optimize_init_params(S,initialvalue)
	p=optimize(S)
	p
	st_matrix("beta_dlwtranslog",p)
	}
end

capture program drop dlw_translog
program dlw_translog, rclass
preserve
**********************************************
sort firm year
mata DLW_TRANSLOG()
end 
use "D:\Nanjing\2019\pollution2\submission\data\boot1",clear
rename frdm firm
gen hangye=int(cic_adj/100)
drop if hangye<13
drop if hangye>42

*keep if hangye==27|hangye==14
capture drop cic2
capture drop cic4 
gen cic2=floor(cic_adj/100)
egen gcic4=group(cic_adj)
egen SDcic=sd(cic2), by(firm)
egen cic_y=group(cic year)
gen Productivity=tfp_acf_new
*----------------------------------------------------------------------
* bootstrap starts from here
*----------------------------------------------------------------------

capture program drop prodboot
program define prodboot, rclass
preserve






* Re-estimate firm-level production function for each sample from bootstrap
do "D:\Nanjing\2019\pollution2\do file\for accre\do file\tfp_combine.do"

xi:heckman ratio1 dyear*, noconstant twostep sel(exp_dummy=  tfp_acf_new clean_dummy_so2_lag clean_dummy_dust_lag clean_dummy_cod_lag dcic2* lnv_lag lnkl_lag dprov* owner* dyear*)  

gen coef1_bs=_b[dyear1]
gen coef2_bs=_b[dyear2]
gen coef3_bs=_b[dyear3]
gen coef4_bs=_b[dyear4]
gen coef5_bs=_b[dyear5]
gen coef6_bs=_b[dyear6]
gen coef7_bs=_b[dyear7]
capture drop exp_dummy2
gen exp_dummy2=.
replace exp_dummy2=exp_dummy*coef1_bs*dyear1 if year==1999
replace exp_dummy2=exp_dummy*coef2_bs*dyear2 if year==2000
replace exp_dummy2=exp_dummy*coef3_bs*dyear3 if year==2001
replace exp_dummy2=exp_dummy*coef4_bs*dyear4 if year==2002
replace exp_dummy2=exp_dummy*coef5_bs*dyear5 if year==2003
replace exp_dummy2=exp_dummy*coef6_bs*dyear6 if year==2004
replace exp_dummy2=exp_dummy*coef7_bs*dyear7 if year==2005
drop  coef1_bs coef2_bs coef3_bs coef4_bs coef5_bs coef6_bs coef7_bs

* For so2 + OLS: column 1 of Table 2
ivreg2 lnso2  (exp_dummy2 so2_new1= treatment_tfp treatment_so2  ventilation)   ///
tfp_acf_new lntariff  i.cic i.prov i.year, partial(i.cic i.prov i.year) cluster (prov_c)

return scalar exp_t2c1=_b[exp_dummy2]
return scalar abt_t2c1=_b[so2_new1]
return scalar tfp_t2c1=_b[tfp_acf_new]
return scalar tar_t2c1=_b[lntariff]


*************************************************************************************************
* For so2+ IV: column 2 of Table 2
ivreg2 lnso2  (exp_dummy2 so2_new1= treatment_tfp treatment_so2  ventilation)   ///
tfp_acf_new lntariff lnkl_lag lnso2_ini i.cic i.prov i.year, partial(i.cic i.prov i.year) cluster (prov_c)

return scalar exp_t2c2=_b[exp_dummy2]
return scalar abt_t2c2=_b[so2_new1]
return scalar tfp_t2c2=_b[tfp_acf_new]
return scalar tar_t2c2=_b[lntariff]
return scalar klr_t2c2=_b[lnkl_lag]
return scalar iei_t2c2=_b[lnso2_ini]



ivreg2 lndust  (exp_dummy2 dust_new1= treatment_tfp treatment_dust  ventilation) ///
  tfp_acf_new lntariff  i.cic i.prov i.year, partial(i.cic i.prov i.year ) cluster (prov_c)

return scalar exp_t2c3=_b[exp_dummy2]
return scalar abt_t2c3=_b[dust_new1]
return scalar tfp_t2c3=_b[tfp_acf_new]
return scalar tar_t2c3=_b[lntariff]

************************************************************************************************************************************
* For dust +IV: column 4 of Table 2
ivreg2 lndust  (exp_dummy2 dust_new1= treatment_tfp treatment_dust  ventilation) ///
  tfp_acf_new lntariff lnkl_lag lndust_ini i.cic i.prov i.year, partial(i.cic i.prov i.year ) cluster (prov_c)

return scalar exp_t2c4=_b[exp_dummy2]
return scalar abt_t2c4=_b[dust_new1]
return scalar tfp_t2c4=_b[tfp_acf_new]
return scalar tar_t2c4=_b[lntariff]
return scalar klr_t2c4=_b[lnkl_lag]
return scalar iei_t2c4=_b[lndust_ini]
***************************************************************************************************************************************


end

capture drop indc
egen indc=group(firm)
drop if indc==.


*------------------------------------------------------------------------
* the following codes generate columns 1, 2, 5, 6  in Table 3
*------------------------------------------------------------------------
bootstrap  exp_t2c1=r(exp_t2c1) abt_t2c1=r(abt_t2c1) tfp_t2c1=r(tfp_t2c1) tar_t2c1=r(tar_t2c1)  ///
 exp_t2c2=r(exp_t2c2) abt_t2c2=r(abt_t2c2) tfp_t2c2=r(tfp_t2c2) tar_t2c2=r(tar_t2c2)  klr_t2c2=r(klr_t2c2) iei_t2c2=r(iei_t2c2)  ///
 exp_t2c3=r(exp_t2c3) abt_t2c3=r(abt_t2c3) tfp_t2c3=r(tfp_t2c3) tar_t2c3=r(tar_t2c3)   ///
 exp_t2c4=r(exp_t2c4) abt_t2c4=r(abt_t2c4) tfp_t2c4=r(tfp_t2c4) tar_t2c4=r(tar_t2c4)  klr_t2c4=r(klr_t2c4) iei_t2c4=r(iei_t2c4),  ///
 reps(50)  seed(1) nodrop cluster(indc) idcluster(newid): prodboot
 
 
 
 
 

 
 
 
 

******************************************************************************************** 
 
 
 
 
 
*-----------------------------------------------------------------------------------------------------
*Part II: The following code generate "bootstrap confidence interval" in Table 3 column 3, 4, 7, 8
*------------------------------------------------------------------------------------------------ ---
   mata:

void GMM_DLW_CD(todo,betas,crit,g,H)
  {
	PHI=st_data(.,("phi"))
	PHI_LAG=st_data(.,("phi_lag"))
	Z=st_data(.,("const","m_lag","l_lag","y1","y2","y3","k","k_lag","exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))
	X=st_data(.,("const","m","l","x1","x2","x3","k","exp_dummy99","exp_dummy00","exp_dummy01","exp_dummy02","exp_dummy03","exp_dummy04","exp_dummy05"))
	X_lag=st_data(.,("const","m_lag","l_lag","x1_lag","x2_lag","x3_lag","k_lag","exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))
	Y=st_data(.,("q"))
	C=st_data(.,("const"))
	

	EXPORT=st_data(.,("exp_dummy99_lag","exp_dummy00_lag","exp_dummy01_lag","exp_dummy02_lag","exp_dummy03_lag","exp_dummy04_lag","exp_dummy05_lag"))

	OMEGA=PHI-X*betas'
	OMEGA_lag=PHI_LAG-X_lag*betas'
	OMEGA_lag_pol=(C,EXPORT,OMEGA_lag)
	g_b = invsym(OMEGA_lag_pol'OMEGA_lag_pol)*OMEGA_lag_pol'OMEGA 
	XI=OMEGA-OMEGA_lag_pol*g_b /*XI=epsilon: recover epsilon_it in equation (23): shock in productivity;*/
	crit=(Z'XI)'(Z'XI)
  }
///
///
void DLW_TRANSLOG()
	{
	initialvalue=st_data(1,("initialConst","initialm","initiall","initialx1","initialx2","initialx3","initialk","initialexp_dummy99","initialexp_dummy00","initialexp_dummy01","initialexp_dummy02","initialexp_dummy03","initialexp_dummy04","initialexp_dummy05"))
	S=optimize_init()
	optimize_init_evaluator(S, &GMM_DLW_CD())
	optimize_init_evaluatortype(S,"d0")
	optimize_init_technique(S, "nm")
	optimize_init_nmsimplexdeltas(S, 0.1)
	optimize_init_which(S,"min")
	optimize_init_params(S,initialvalue)
	p=optimize(S)
	p
	st_matrix("beta_dlwtranslog",p)
	}
end

capture program drop dlw_translog
program dlw_translog, rclass
preserve
**********************************************
sort firm year
mata DLW_TRANSLOG()
end 

use "D:\Nanjing\2019\pollution2\submission\data\boot_app2",clear
*rename party_id frdm

*
*----------------------------------------------------------------------
* bootstrap starts from here
*----------------------------------------------------------------------
*keep if _merge==2|_merge==3
capture program drop prodboot
program define prodboot, rclass
preserve






* Re-estimate firm-level production function for each sample from bootstrap



xi:heckman ratio1 dyear*, noconstant twostep sel(exp_dummy=  tfp_acf_new clean_dummy_so2_lag clean_dummy_dust_lag clean_dummy_cod_lag dcic2* lnv_lag lnkl_lag dprov* owner* dyear*)  

gen coef1_bs=_b[dyear1]
gen coef2_bs=_b[dyear2]
gen coef3_bs=_b[dyear3]
gen coef4_bs=_b[dyear4]
gen coef5_bs=_b[dyear5]
gen coef6_bs=_b[dyear6]
gen coef7_bs=_b[dyear7]
capture drop exp_dummy2
gen exp_dummy2=.
replace exp_dummy2=exp_dummy*coef1_bs*dyear1 if year==1999
replace exp_dummy2=exp_dummy*coef2_bs*dyear2 if year==2000
replace exp_dummy2=exp_dummy*coef3_bs*dyear3 if year==2001
replace exp_dummy2=exp_dummy*coef4_bs*dyear4 if year==2002
replace exp_dummy2=exp_dummy*coef5_bs*dyear5 if year==2003
replace exp_dummy2=exp_dummy*coef6_bs*dyear6 if year==2004
replace exp_dummy2=exp_dummy*coef7_bs*dyear7 if year==2005
drop  coef1_bs coef2_bs coef3_bs coef4_bs coef5_bs coef6_bs coef7_bs







*************************************************************************************************

 ivreghdfe lnso2  (Export Abatement_so2= tariff_prod ventilation_prod  ventilation) Productivity Tariff lnKL , a(indc_firm year)

return scalar exp_t1c2=_b[Export]
return scalar abt_t1c2=_b[Abatement_so2]
return scalar tfp_t1c2=_b[Productivity]
return scalar tar_t1c2=_b[Tariff]
return scalar klr_t1c2=_b[lnKL]


capture drop ub
egen ub=pctile(lndust), p(97)
 ivreghdfe lndust  (Export Abatement_dust= tariff_prod ventilation_prod  ventilation) Productivity Tariff lnKL if lndust<ub, a(indc_firm year)

return scalar exp_t2c2=_b[Export]
return scalar abt_t2c2=_b[Abatement_dust]
return scalar tfp_t2c2=_b[Productivity]
return scalar tar_t2c2=_b[Tariff]
return scalar klr_t2c2=_b[lnKL]
********************************************************************************************************


*----------------------------------------------------------------------------------------------------
* Column 4 of Table2
ivreghdfe lnso2_d  (exp_dummy2 so2_new1_d= treatment_tfp treatment_so2  ventilation)  ///
 tfp_acf_new_d lntariff_d lnkl_d, a(cic prov year ownership lndust_ini)


return scalar exp_t3c2=_b[exp_dummy2]
return scalar abt_t3c2=_b[so2_new1_d]
return scalar tfp_t3c2=_b[tfp_acf_new_d]
return scalar tar_t3c2=_b[lntariff_d]
return scalar klr_t3c2=_b[lnkl_d]

*----------------------------------------------------------------------------------------------------
* Column 8 of Column 2
ivreghdfe lndust_d  (exp_dummy2 dust_new1_d= treatment_tfp treatment_dust  ventilation)  ///
 tfp_acf_new_d lntariff_d lnkl_d , a(cic prov year ownership lndust_ini) 
return scalar exp_t4c2=_b[exp_dummy2]
return scalar abt_t4c2=_b[dust_new1_d]
return scalar tfp_t4c2=_b[tfp_acf_new_d]
return scalar tar_t4c2=_b[lntariff_d]
return scalar klr_t4c2=_b[lnkl_d]


end

capture drop indc2
egen indc2=group(firm year)
drop if indc2==.


**************************************************************************************************
* The following code generate columns 3, 4 ,7, 8 in Table 3.
*************************************************************************************************
bootstrap   exp_t1c2=r(exp_t1c2) abt_t1c2=r(abt_t1c2) tfp_t1c2=r(tfp_t1c2) tar_t1c2=r(tar_t1c2)  klr_t1c2=r(klr_t1c2) ///
exp_t2c2=r(exp_t2c2) abt_t2c2=r(abt_t2c2) tfp_t2c2=r(tfp_t2c2) tar_t2c2=r(tar_t2c2)  klr_t2c2=r(klr_t2c2) ///
exp_t3c2=r(exp_t3c2) abt_t3c2=r(abt_t3c2) tfp_t3c2=r(tfp_t3c2) tar_t3c2=r(tar_t3c2)  klr_t3c2=r(klr_t3c2) ///
exp_t4c4=r(exp_t4c2) abt_t4c2=r(abt_t4c2) tfp_t4c2=r(tfp_t4c2) tar_t4c4=r(tar_t4c2)  klr_t4c4=r(klr_t4c2), ///
 rep(50)  seed(5) nodrop : prodboot
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 